************************************************
* Sleep Project - Pedro Bessone, Gautam Rao, Heather Schofield, Frank Schilbach, and Mattie Toma
* Purpose: Replicates main paper, Figure 5 (Impacts on Nighttime and Nap Sleep)
* Last edited: 07 May 2021
************************************************


	************************************************	
	*0. Initial Setup
	************************************************

	clear all
	set more off
	set matsize 800
	
************************************************	
*1. First stage: 
************************************************

	*Prep data
	
		*Analysis base file
		use "$d/firststage_dataset.dta", clear
		
		gen treat_int = treat_pool*treat_nap
		gen treat_int_s = treat_s*treat_nap
		gen treat_int_s_i = treat_s_i*treat_nap
		
		gen treat_pool_cell = treat_pool*(1-treat_nap)
		gen treat_s_cell = treat_s*(1-treat_nap)
		gen treat_s_i_cell = treat_s_i*(1-treat_nap)
		gen treat_nap_cell = treat_nap*(1-treat_pool)
		gen treat_int_cell = treat_int 
		gen treat_int_s_cell = treat_int_s
		gen treat_int_s_i_cell = treat_int_s_i
		
		replace nap_time_mins = .1 if pid == 5130 & day_in_study == 27 //raw actigraph data incorrectly recorded as half hour nap
	
		* Mean Daily Graphs
		replace  nap_time_mins = 0 if post_treatment == 0

	* Mean Daily Graphs - Fully Disaggregated

		gen treatment_disag = 1 if treat_s_cell == 1
		replace treatment_disag = 2 if treat_s_i_cell == 1
		replace treatment_disag = 3 if treat_nap_cell == 1
		replace treatment_disag = 4 if treat_int_s_cell == 1
		replace treatment_disag = 5 if treat_int_s_i_cell == 1
		replace treatment_disag = 0 if mi(treatment_disag)

	
	* Mean Daily Graphs - Participants in one cell only


		local x=0
		
		foreach var in Sleep_Night nap_time_mins act_inbed epoch_sleepefficiency {
					
			preserve
					
			drop if day_in_study == 1					
			local max_day = 28
					
			if "`var'" == "nap_time_mins" {
				drop if day_in_study == 28
				local max_day = 27
			}
								
		*Sample sizes
		quietly summarize `var'
		local ss = r(N)

		*Number of PIDs
		distinct pid
		local pids = r(ndistinct)
				
		*Gen treatment and control mean by day
							
				*Gen normalized means
				bysort day_in_study treatment_disag: egen mean_`var'_tc = mean(`var')

		*Line graph of treatment group means	
		if "`var'" == "nap_time_mins" {
		
			#d ;
				twoway	
				(line mean_`var'_tc day_in_study if treatment_disag == 0, lcolor(black) lpattern(dash) lwidth(medthick) sort)
				(line mean_`var'_tc day_in_study if treatment_disag == 3, lcolor(green) lpattern(longdash) lwidth(medthick) sort),
				ytitle("`t`x''")
				xtitle(Day in Study)
				xlabel(2(2)`max_day', labsize(small)) ysc(r(0 16)) ylab(0(4)16, labsize(small) gmin gmax)
				legend(order(1 "Control" 2 "Nap Only" ) symxsize(*1) size(*1) region(col(white)) col(3))
				xline(8.5, lcolor(black))
				graphregion(fcolor(white))
				name(graph`x', replace)
				
			;
			#d cr
		
		}
		else if "`var'" == "epoch_sleepefficiency" {
		
			#d ;
				twoway	
				(line mean_`var'_tc day_in_study if treatment_disag == 0, lcolor(black) lpattern(dash) lwidth(medthick) sort)
				(line mean_`var'_tc day_in_study if treatment_disag == 1, lcolor(red) lpattern(longdash) lwidth(medthick) sort)
				(line mean_`var'_tc day_in_study if treatment_disag == 2, lcolor(blue) lpattern(shortdash) lwidth(medthick) sort),
				ytitle("`t`x''")
				xtitle(Day in Study)
				xlabel(2(2)`max_day', labsize(small)) ysc(r(0.6 0.8)) ylab(0.6(0.05)0.8, labsize(small) gmin gmax)
				legend(order(1 "Control" 2 "Encouragement Only" 3 "Incentives Only") symxsize(*1) size(*1) region(col(white)) col(3))
				xline(8.5, lcolor(black))
				graphregion(fcolor(white))
				name(graph`x', replace)
				
			;
			#d cr
		
		}
		
		else {
		
		
			#d ;
				twoway	
				(line mean_`var'_tc day_in_study if treatment_disag == 0, lcolor(black) lpattern(dash) lwidth(medthick) sort)
				(line mean_`var'_tc day_in_study if treatment_disag == 1, lcolor(red) lpattern(longdash) lwidth(medthick) sort)
				(line mean_`var'_tc day_in_study if treatment_disag == 2, lcolor(blue) lpattern(shortdash) lwidth(medthick) sort),
				ytitle("`t`x''")
				xtitle(Day in Study)
				xlabel(2(2)`max_day', labsize(small))
				ylabel(`yl`x'', labsize(small) gmin gmax)
				legend(order(1 "Control" 2 "Encouragement Only" 3 "Incentives Only") symxsize(*1) size(*1) region(col(white)) col(3))
				xline(8.5, lcolor(black))
				graphregion(fcolor(white))
				name(graph`x', replace)
				
			;
			#d cr		
		
		
		}
				
				
				
			*Save graphs

				graph export "$om/Figures/Figure5/Mean_Daily_`var'_only_one_cell.eps", replace	 
			
				
			
		local x = `x' + 1	
				
		restore
	}	
	
